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ABSTRACT 

We consider the pulsational properties of white dwarf star models with temperatures appropriate for the ZZ Ceti 
instability strip and with masses large enough that they should be substantially crystallized. Our work is motivated 
by the existence of a potentially crystallized DAV, BPM 37093, and the expectation that digital surveys in progress 
will yield many more such massive pulsators. 

A crystallized core makes possible a new class of oscillations, the torsional modes, although we expect these 
modes to couple at most weakly to any motions in the fluid and therefore to remain unobservable. The /9-modes 
should be affected at the level of a few percent in period, but are unlikely to be present with observable amplitudes 
in crystallizing white dwarfs any more than they are in the other ZZ Ceti's. Most relevant to the observed light 
variations in white dwarfs are the g-modes. We find that the kinetic energy of these modes is effectively excluded 
from the crystallized cores of our models. As increasing crystallization pushes these modes farther out from the 
center, the mean period spacing (AP) between radial overtones increases substantially with the crystallized mass 
fraction, Mc, /M*. In addition, the degree and structure of mode trapping is affected. The fact that some periods 
are strongly affected by changes in the crystallized mass fraction while others are not suggests that we may be 
able to disentangle the effects of crystallization from those due to different surface layer masses. 

Subject headings: dense matter — stars: oscillations, evolution — white dwarfs 



1. ASTROPHYSICAL CONTEXT 

The theoretical study of pulsating crystalline objects extends 
ma ny years into the past. One of t he first numerical studies was 
by Sherman, Jarosch, & Pekeris (1959), who modeled global 



oscillations of the Earth. Their main interest was in fitting the 
oscillation period of 57 minutes which was excited by the Kam- 
chatka earthquake of 1952. In the process, they examined how 
the central density in their models allowed them to match the 
periods of other oscillation modes which were also observed to 
be excited by the earthquake. 
In an astrophysical context. 



Hansen & Van Horn (1979) 



treated oscillations in white dwarf models with a crystalline in- 
ner core. Since it was known that 1 Mq models with Tf^g ~ 
10,0 00 K were in the process of crystallizing (Lamb & Van 



Horn 1975; [Van Horn & Savedoff] 1976), [Hansen & Van Horr^ 
self-consistently treated the response of the crystalline core to 
the pulsations. Their main interest was in explaining the ob- 
served period ranges of the ZZ Ceti's in terms of low-radial 
order oscillations. They found that the ^-mode periods were 
decreased by the presence of crystallization, contrary to our 

present findings. 

McDermott, Van Horn, & Hansen (1988) treated oscillations 
in neutron star models with a fluid core, a solid crust, and a 
thin surface fluid "ocean." They considered neutron star oscil- 
lations as a possible explanation for the observed irregularities 
in the timing of subpulses from radio pulsars, and as a source 
of the observed periodicities in many X-ray burst sources. They 
found g-modes which were trapped in the cores of their models, 
as well as t hose which were tra pped in the surface oceans. 
Finally, Bildsten & Cutler (1995) considered g-mode os- 



tems. They found a good match to these frequencies for low 
order, £ = I g-modes. 

Why, then, does this problem need to be re-examined in the 
context of white dwarf stars? As is often the case, new observa- 
tions and new circumstances have again made this problem one 
worth considering, but in more detail than the general analyses 
of the past. For example, the pioneering calculations of Hansen 
& Van Horn (1979) were focussed primarily on the range of 
normal mode periods which are possible given a crystallized 
core, not with the details of how the periods of high overtone 
g-modes are affected at the level of 5-10%. At the time, there 
were no known high-mass white dwarf pulsators, and precise 
mode identifications for any pulsating white dwarf had yet to 
be attempted. 

That situa tion changed with the discovery of pulsations in 
BPM 37093 dKanaanetai] 1992), a high-mass ZZ Ceti star (see 
Figu re [l|) which should be substantially crystallized ( Winget 
et al. 1997); depending on the C/O ratio in its core, it should 
be between 50% and 90% crystallized by mass. Depending on 
the details of its nuclear history, its core could be composed 
of even heavier elements such as Ne ( Iben 1991), which would 
imply that it is more than 90% crystallized (Winget et al. 1997). 



The Whole Earth Telescope (WET) examined this target in the 
Spring of 1998 and found at least 8 independ ent freq uencies, 
three of which had been previously seen by |Kanaan| (1996). 



dilations in the thin surface oceans of accreting neutron star 
models. Their aim was to explain the observed 5-7 Hz quasi- 
periodic oscillations in the brightest accreting neutron star sys- 



Thus, the potential to perform asteroseismology on this object 
requires us to make a more detailed theoretical investigation of 
the properties of crystallized pulsators. 

One hope is that we will be able to independently deter- 
mine the crystallized mass fraction Md/M*, and thereby pro- 
vide a direct test of the theor y of c rystallizati on, no w nearly 
four decades old (Abrikosov 1960; Kirzhnits] 1960; Salpetei 
1961). This subject is relevant to the astronomical commu- 
nity at large, since phase separation of C and O during crystal- 
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Fig. 1 . — The position of BPM 37093 relative to the other ZZ Ceti's in Bergeron et al. 1995 as a function of T^ff and M^/Mq . The Hnes coiTespond to constant 
amounts of crystallized mass fraction assuming a pure oxygen core. If BPM 37093 has an oxygen core it should be ~ 90% crystallized, and for a carbon core, ~ 
50% crystallized. 



lization, and, indeed, crystallization itself, represent the largest 
sources of systematic uncertainties in the age of the local Galac- 
tic disk as derived from the white dwarf luminosity function. In 
addition, understanding the internal structure of white dwarfs 
may prove vi tal in fitting cosm ological models to Supernova la 



(SNIa) data (3arnavich et al. 1998), so that systematic differ- 
ences in the absolute magnitudes of the SNla may be corrected 
for the evolutionary diff erences in the SN progenitors (H oflich, 
Wheeler, & Thielemann 1998). 

Finally, digital surveys now in progress promise to add con- 
siderably to the presently known number of cool white dwarfs. 
For instance, the Sloan Digital Sky Survey (3unn 1995) should 
increase this number by a factor of approximately 20, with the 
result that we may have 20 such stars with which to test the 
theory of crystallization. 

2. EVOLUTIONARY MODELS 

The basis for our equilibrium models is an updated version 
of WDEC, the White Dwarf Ev olution ary Code, as described 



in Lamb & Van Horn 



(1975) and Wood (1990, 1992). Here we 



present only a brief summary of the input physics in our mod- 
els, with references provided for a more complete description. 

In the cores of our models we use the Lamb equation of state 
(EOS) dLamt 1974) and in the envelope s we use the tabular 
EOS of [Fontaine, Graboske, & Van Horn (1977). We employ 
the additive volume technique to interpolate between pure com- 
positions for the Carbon/Oxygen mixture in the core and the 
Hydrogen/Helium/Carbon mixture in the envelope. The chem- 
ical profiles of the composition transition zones in the enve- 
lope are treated with an adaptation of the method of A rcoragi & 
Fontaine (1980). Essentially, these profiles mimic those which 
would be obtained in diffusive equilibrium, but contain addi- 
tional parameters controlHng the thickness of the transition re- 



gions ( [Bradley, Winget, & Wood) 1993). 

The question of crystallization and our treatment of it is cen- 
tral to our analysis. For a model with a given mass, T^s, and 
composition, the Lamb EOS does return a unique answer for 
the degree of crystallization: the critical value of F is given 
by Fcr ~ 160, where F = Z^e^/ {r)kBT is the ratio of Coulomb 
energy between neighboring ions to each ion's kinetic energy. 
More recent calculati ons indicate a some what higher value for 
this ratio, Fcr — 180 ( Dgata & Ichimari^ 1987). Our approach 
is to compute equilibrium models using WDEC and the Lamb 
EOS, which results in models with a self-consistently computed 
value of the crystallized mass fraction. When we perform a pul- 
sational analysis of these models, however, we take the crystal- 
lized mass fraction to be a free parameter, in hopes of using 
asteroseismology to place constraints on the degree of crystal- 
lization. The underlying assumption here is that two equilib- 
rium models which differ only in the degree of crystallization 
have virtually identical pressure, density, and temperature pro- 
files. 

While this is not a physically self-consistent procedure, it 
is justifiable for two reasons. First, the main physical effect 
of crystallization is the release of latent heat; this provides the 
models with an additional energy source, which means that at 
a given Teff they are older. This clearly has no effect on the 
pulsational properties, which depend only on the structural pa- 
rameters of a given equilibrium model. Second, the density 
change at crystallization is quite small, Sp/p ^ 10~^ ( Lamb & 
Van Horn 1975), so the difference in, for example, p(r), P{r), 
and T{r) between two models which differ only in the amount 
of crystallization is accordingly quite small. As we show in 
section 5.2, the effect of crystallization upon g-mode pulsations 
can be accurately taken into account through a modified bound- 
ary condition at the assumed solid/fluid interface. 
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3. ASYMPTOTIC NONRADIAL OSCILLATION THEORY 

Stars which are fluid (uncrystalHzed) can undergo nonradial 
motions which have been labeUed g- and p-mode osciflations. 
In the linear limit, these modes of oscillation are spheroidal, 
with the Eulerian perturbations of variables such as the density 
and pressure having the angular spatial dependence of a single 
spherical harmonic (i.e., p', p' cx Y'^'iO, </>)). From a local anal- 
ysis, the radial wavenumber A:,, is given by 



modes, in this case), and (2) the pre-existing p- and ^-modes 
are modified. We now treat these cases separately. 

4. 1 . The Torsional Modes 

The torsional modes, or f-modes, are very special nonradial 
modes characterized by zero radial displacement and zero com- 
pression, i.e., and V • ^ both vanish, where ^ is the displace- 
ment vector. The dispersion relation for these modes is 



1 



(1) 



where a is the angular frequency of the mode, c, is the sound 
speed, L\ = £(£+ l)c^/r^ is the square of the Lamb/acoustic 
frequency, r is the radi al variable, a nd A^^ is the famed Brunt- 
Vaisala frequency (see Unno et al. 1989 for a more complete 
discussion). From the above formula, we see that a mode is 
propagating in a region (i.e., has kf. > 0) if cr^ > L\,N'^, or 
< LjjN^. Thus, the modes separate cleanly into two classes: 



(2) 



p-modes 



where = /i/p is the square of the shear velocity, p is the shear 
modulus, and T/^ = [£(£+ l)-2]v^ is the "torsional frequency." 
They propagate in the region defined by > Tg, and their fre- 
quency spectrum is equally spaced, as is the case with /?-modes, 
with 



> Lj,N^, "high-frequency limit" 



displacements become vertical near the surface 



^-modes 



"low-frequency limit" 

^:lNdr/rY 



Jdr/v,' 

As we might expect, the k=l period for these modes goes like 
Ri,/vs, the crossing time for a shear wave. 

In Figure | , we show a propagation diagram for f-modes with 
£ = 2, using the same white dwarf model as in Figure ^. If we 
imagine a model which is 90% crystallized, then the f-mode 
can potentially propagate anywhere inside the 90% mass point 
in the model. If the mode is an £ = 2 mode, then its region 
of propagation is restricted further to the region for which its 
frequency is greater than the torsional frequency, i.e., cr^ > T2- 
For a 1 sec mode, this corresponds to the part of the horizontal 
dotted line which lies to the left of the 90% point in Figure |[ 
We note that for all the £ = I modes, we have T^ = 0, so these 
modes propagate throughout the entire crystallized region. 

Observable consequences of the f-modes, if any exist, are 
difficult to identify. The longest period f-modes should have 
periods ^ 20 sec, which is too short to explain the observed 
oscillations in the ZZ Ceti's of lOO's of seconds.' In addition, 
these modes should not be able to couple (in the linear limit) to 
the fluid at the solid/fluid interface, so these oscillations should 
be unable to propagate from the crystalline core through the 
fluid to the surface. Also, the different angular structure of 
the f-modes should make any nonlinear coupling between these 
modes and the ordinary p- and g-modes very weak; to the first 
nonlinear order this coupling will be zero. As a result, we ex- 
pect these modes to be unobservable unless crystallization has 
proceeded out into the photosphere. The oldest known white 
dwarfs in the Galaxy are not yet cool enough for this to have 
occurred. We therefore turn our attention to the p- and g-modes. 

4.2. The Spheroidal Modes 

4.2.1. p-modes 

For pressure waves traveling in a solid medium, the velocity 
Vp is given by 

, A + 2/i 

P 

where A = Fif -2/i/3, p and P are the density and pressure, 
respectively, and Fi is the usual adiabatic exponent (Landau & 

'Such short timescale oscillations of lO's of seconds may be relevant for accreting white dwarfs in cataclysmic variable systems, although the heating due to 
accretion may preclude the presence of substantial crystallization in these objects. 



displacements become horizontal near the surface 

Here rl and r2 are the inner and outer classical turning points, 
respectively, at which = for a given a. We see that in 
the asymptotic limit the /?-modes are uniformly spaced in fre- 
quency as a function of radial order k, while the ^-modes are 
uniformly spaced in period. 

A useful diagnostic for the frequency spectrum of a given 
white dwarf model is the propagation diagram, an example of 
which is shown in Figure |[ where we have labelled the high 
and low frequency domains of the p- and g-modes; the model 
is of a 1.1 Mq white dwarf with T^ff = 12,200 K. The horizon- 
tal axis is given in terms of the radial variable ln(r//?), which 
is the natural logarithm of the radius divided by the pressure, 
where both r and p are given in cgs units. This radial variable 
has the desirable property that it increases monotonically out- 
ward from the center and provides increased resolution in both 
the core (where r approaches zero) and in the envelope (where 
p approaches zero). Along the upper axis we display the more 
commonly used radial variable -log(l -Mr/M*), which may be 
more easily related to the structural parameters of the models. 

We note that the bumps in A^^ and L| correspond to the 
C/O, He/C, and H/He transition zones. For instance, using 
the upper axis to obtain estimates of -log(l -M^/M^^), we 
see that the outer two transition zones have Mu/Mj, ^ 10"^, 
Mee/A^i ~ 10"^, which are in fact the values assumed in these 
models. 

4. THE EFFECT OF A CRYSTALLINE CORE 

How does a crystalline core affect the oscillations of a star? 
As we introduce a solid core into our models, two things oc- 
cur: (1) A new class of modes appears (the torsional/toroidal 
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Fig. 2. — A propagation diagram showing Af^ and Lj as a function of In (r/p) (lower axis) and -log(l -Mr /Mi,) (upper axis); the center is on the left and the 
surface is on the right. The region of propagation of a 600 second g-mode is shown. The vertical dashed lines are labelled by the percent mass which is interior to 
these regions, i.e., the 90% line indicates the boundary at which 90% of the mass of the model is inside this point. We see that a model which is this crystallized 
now has an inner turning point for g-mode propagation considerably farther out than in the uncrystallized case. 




Fig. 3. — A propagation diagram for t-modes for i = 2. The t-modes propagate only in the crystallized region, e.g., only to the left of the 90% crystalUzed line for 
a 90% crystalUzed model. We note that for ^ = 1 we have = 0, so the modes propagate throughout the entire crystallized region in this case. 
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Lifshitz 1975). If we treat the non-zero as a perturbation, we 
find that 



Sv„ 



2/i 

3riP' 



where Svp is the change in Vp due to the finite shear modulus. In 
the cores of our 1.1 Mq models, we typically find fJ./p ^ 0.01. 
Thus, p-mode periods are affected at the level of only a few per- 
cent by the presence of a crystalline lattice. They are therefore 
of no more interest than are ordinary /9-modes in the context of 
the observed pulsations of the DAV and DBV white dwarfs. 

4.2.2. g-modes 

We concentrate the remainder of our analysis on the g- 
modes, since these are the modes which are believed to be re- 
sponsible for the observed pulsations in the white dwarf vari- 
ables. Because ^-modes have large shears associated with their 
fluid motions, we expect the nonzero shear modulus of the 
solid to have a significant effect on them. Qualitatively, we 
may ask when the return force due to a finite shear modulus is 
approximately equal to the return force normally experienced 
by flu id elements in the absence of such shear (e.g., B ildsten & 
Cutler 1995). Algebraically, the shear return force is equal to 
or exceeds the ordinary return force of the fluid when 



> 1, 



where h = P/\dP/dr\ is a pressure scale height. 



In our models, we find that 



> 10 , which indicates 



that the g-modes are completely altered in the crystallized re- 
gion. Thus, a g-mode which is propagating in the fluid region 
will find a complete mismatch as it attempts to propagate into 
the crystallized region. We therefore expect nearly complete 
reflection of the g-mode at such a boundary, with the result that 
the g-modes are essentially confined to the fluid regions of our 
models. 

5. NUMERICAL ANALYSIS 

5.1. The Global Solution 

We now examine the above assertion and offer a numerical 
justification for it. Our approach is based on the work of Hansen 
& Van Horn (1979); we treat the "global" problem in that we 
allow the solid cores of our models to respond to the oscilla- 
tions. We have used the Cowling approxi mation to simplify the 
pulsation equations, as was also done in ^ansen & Van Horq 
(1979). Since g-modes in white dwarfs are primarily envelope 
modes, this is an excellent approximation and hardly affects 
the accuracy of our calculated periods; even k= I, £ = 2 modes 
have pe riods which are only affected at the level of 0.2% (M ont- 
gomery 1998). The details of the rest of the global treatment are 
summarized in Appendix ^ where we describe the oscillation 
variables, the equations which they obey, the central boundary 
conditions, and the connecting conditions at the solid/fluid in- 
terface. 

In Figure Q we plot the radial and horizontal displacements 
of a 378.4 sec, £ = 1 mode; the model is again that of a 1.1 Mq 
white dwarf with Teff = 12,200 K, which is assumed to be 50% 
crystallized. As is true of all the g-modes we have examined, 
the amplitude of the fluid motions is decreased by 3 orders 
of magnitude in the solid as compared to the fluid. One other 
feature of the oscillations is that the horizontal displacement is 
discontinuous at the solid/hquid interface. In the approximation 



of zero viscosity and laminar flow, the fluid is free to slide over 
the solid surface. In reality, a turbulent boundary layer would 
probably form in this region, which would tend to dissipate the 
pulsation energy. 

The kinetic energy density depends on the square of the dis- 
placement, so it is attenuated by ~ 6 orders of magnitude in the 
solid core. Since the kinetic energy is an indicator of how a 
given mode samples the different regions of a model, we con- 
clude that it is a very good approximation to treat the g-modes 
as excluded from the solid cores of our models. In the fol- 
lowing section we will demonstrate the validity and the self- 
consistency of this approach. 

5.2. The "Hard-Sphere" Boundary Condition 

As suggested in the previous section, we may be able to 
reproduce the effects of crystallization on g-mode pulsations 
merely by applying a hard-sphere boundary condition at the 
solid/liquid interface. By this we mean that the radial dis- 
placement is set to zero (^,- = 0) and the horizontal displace- 
ment is left to be arbitrary. In addition, the boundary condition 
on the gravitational potential and its derivative are the same as 
for the uncrystallized case, as we show in appendix ^ these 
"hard-sphere" calculations are not in Cowling approximation 
and therefore solve the full fourth-order adiabatic equations in 
the fluid region. Using the hard-sphere boundary condition has 
the advantage that the resulting problem is much easier to treat, 
both in terms of speed and convergence. 

We have calculated the fractional difference between peri- 
ods calculated with the "hard-sphere" approximation and those 
calculated with the "global" treatment. Using a fiducial model 
with = 1.1M0, Teff = 12,200 K, and assuming 90% crys- 
tallization by mass, we have examined all ^ = 1 and 2 periods 
between 50 and 1000 sec. We find that the fractional difference 
in periods is less than 1 part in 10"^, and that the absolute error 
in the calculated periods never exceeds 0.05 sec. We there- 
fore conclude that the "hard sphere" boundary condition at the 
solid/fluid interface accurately represents the physics of g-mode 
os cillations in models with crystalline cores. B ildsten & Cut- 
ler (1995) found exactly the same approximation to be valid in 
their treatment of g-modes in the surface oceans of accreting 
neutron star models. 

Before proceeding to the detailed numerical calculations, we 
wish to convince the reader that crystallization will have a mea- 
surable effect on the periods. In Figure |] we have plotted the 
kinetic energy per unit x = ln{r/ p), so that the area underneath 
the curve represents the weight of each region's contribution to 
the total kinetic energy as a function of x. The vertical dashed 
lines indicate different mass points in this model. For instance, 
if the model is 90% crystallized, then the kinetic energy to the 
left of the 90% line is eliminated from the mode. By visual in- 
spection, this is of order 10% of the kinetic energy in the mode, 
so we might well expect that the period of this mode is affected 
at the 10% level. In fact, we will see in the next section that the 
periods can be shifted by even larger amounts. 

6. THE G-MODE PERIODS AS A FUNCTION OF Mcr/M* 
6.1. Asymptotic Relations 



The kinetic energy argument in section 5.2 leads us to expect 
that the g-mode periods will change measurably as the crys- 
tallized mass-fraction increases from to 90%. With this in 
mind, we re-examine the asymptotic formulae for g-mode peri- 
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Fig. 4. — The log of the absolute values of the radial (upper panel) and horizontal (lower panel) displacements as a function of ln(r/p). Note that is continuous 
at the solid/crystal interface at In (r/p) ~ -36.7, but that is not. The magnitudes of both and are reduced by ~ 3 orders of magnitude as they penetrate the 
solid region. 




Fig. 5. — The kinetic energy per unit x = ln(r/p). The vertical dashed lines are labelled with the mass-fraction of the model interior to the given point. The 
mode shown is an £ = l,i = 25 mode with a period of 673.4 sec. The equilibrium model is a 1.1 Mq model with T^ff at 12,200 K, and the surface layer masses are 
MHe/M* = 10"' andMn/M, = lO"'. 
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ods (e.g., LJnno et al. 1989) 



{AP): 



VW+T) 



Ndr/r 



(3) 



where we have written (AP) for the mean period spacing be- 
tween consecutive radial orders. Since the g-modes are ex- 
cluded from the crystallized region, the inner turning point rl 
is now a function of M^^/Mi,. As we allow the model to crys- 
tallize while holding all other structural parameters constant, rl 
moves outward, so the integral in equation]^ decreases, with the 
result that {AP) and Pk both increase. 

As a heuristic tool, we would like to plot the "region of pe- 
riod formation," which would tell us visually the weight which 
the different regions of the star have in determining the period 
of a mode. This proble m has been examined several t imes in 
the past, for exa mple by Kawaler, Winget, & Hansen (1985), 
Schwank (1976), Goosens & Smeyers (1974), and originally 
by l^pstein| (1950). 

ib simplify matters, we examine this weight function in the 
asymptotic limit of high k and I. For ^-modes, we have to be 
content to determine a "region of frequency formation." From 
asymptotic theory, we find that the relative contribution to the 
total frequency per unit radius is 

da 

dr r ' 

which depends only on and r. In order to expand the radial 
axis in both the center and the envelope and to make the result- 
ing functions easier to examine, we choose x = ln(r/ p) as our 
radial coordinate. Then the above relation becomes 



da 
dx 



N 

TTy' 



(4) 



where V = Tigr/c^. We emphasize that the appearance of the 
sound speed in the variable V is purely a result of the above 
radial coordinate change, and does not reflect a dependence of 
g-mode frequencies on c,. 

In Figure |[ we plot da/dx versus x for a 1.1 Mq, Teff = 
12,200 K model with Mue/M^ = 10"^ and Mh/M^ = 10"^ The 
three spikes in da/dx correspond to the composition transition 
zones of O/C, C/He, and He/H. From inspection of this figure, 
we would expect the C/He transition zone to have the least ef- 
fect on the g-mode periods, whereas the He/H transition zone 
in the e nvelope should have the largest effect. Numerically, 



Bradley (1993) has found this to be the case, with the period 



spacing and mode trapping being most sensitive to the hydro- 
gen layer mass and least sensitive to the thickness of the helium 
layer Physically, this is due to the fact that the He/H transition 
zone, since it is closest to the surface, is the least degenerate, 
so it has the largest thermal contributions to the Brunt-Vaisala 
frequency. In addition, it is the only zone in which there is a 
contrast in the atomic weight per electron, In going from 
He to H, He goes from 2 to 1; for the O/C and C/He zones, 
lie =2 for both chemical species in the transition zone. 

It is worth commenting on the similarities between the dis- 
tribution of kinetic energy in Figure [s] and the shape of the g- 
mode period formation region in Figure |[ The kinetic energy 
plot is for a numerically calculated £= l,k = 25 mode, whereas 
the period formation region is in the high k hmit. The value of 



N/(l + V) in Figure ^ should correspond to to the wavelength 
of oscillations as a function of x in Figure |] This is in fact the 
case, since we see that peaks in Figure ^ correspond to rapid 
spatial oscillations in the kinetic energy density. Similarly, the 
small value ofA^/(l+y)in Figure ^ for x in the range of - 1 to 
-4 results in a longer spatial wavelength in the oscillations of 
the kinetic energy density in Figure || at this value of x. 

We also note that the overall envelope of the kinetic energy is 
similar in shape to Figure |[ With this in mind, we calculate the 
kinetic energy distribution in the asymptotic limit. The kinetic 
energy dE in a shell dr is given by 

dE w pr^dr[^f + £(£+l)^f/\ 

where we have used the fact that ^ for g-modes. If we 
now substit ute for f/, the asymptotic value for it taken from 
Unno et al. (1989), then we find that the envelope of the ki- 



netic energy density varies like 

dE _ 
dr 



N 



(5) 



Thus, we see that in the asymptotic limit the kinetic energy 
samples the model in the same way as does the frequency for 
g-modes. 

Although the above discussion might give the impression that 
modes of this radial overtone number are safely in the asymp- 
totic limit, such is not the case. For the mode in Figure || it 
may be treated in the asymptotic limit in the region between 
X = -20 and x = -5, i.e., the variations of its amplitude are small 
compared to its wavelength. However, the composition transi- 
tion zone at X ~ -23 provides a much more rapid spatial vari- 
ation than the wavelength of this mode. Depending on the de- 
tails of how the mode interacts with this feature, it will be par- 
tially transmitted and partially reflected at this boundary. Thus, 
the amplitudes of the mode on each side of a transition zone 
will not in general be given by the asymptotic theory. In other 
words, the effect of a transition zone is to enhance the ampli- 
tude of a mode on one side of a transition zone relative to its 
amplitude on the other side. This effect is generically known as 
"mode trapping," although in the context of white dwarfs this 
term usually denotes a mode which has an enhanced amplitude 
in the outer surface layer, i.e., the H layer for DAV's. Return- 
ing to the mode in Figure ^ neighboring modes which differ 
from it by only ± 1 in still have somewhat different distribu- 
tions of kinetic energy between the different transition zones. 
Thus, we cannot consider these modes to be globally described 
by asymptotic theory. 

6.2. Numerical Results 

We now wish to make a comparison between the functional 
form of the period spacing implied by equation || and that de- 
rived from direct numerical calculations. To do this, we nor- 
malize (AP) to the average period spacing in the uncrystal- 
lized case, denoted by (AP)o. Such a comparison is shown 
in Figure where the solid line gives the analytic relation 
and the filled circles are the result of a numerical pulsational 
analysis of £ = 2 periods between 500 and 1000 sec. We have 
made the model, a 1.1 Mq C/O core model with T^s = 12,200 
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Fig. 6. — The frequency (period) formation region for g-modes in a 1.1 Mq model with T^g = 12,200 K. The three spikes are all composition transition zone 
features, which from left to right are due to the O/C, C/He, and He/H transition zones. 
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Fig. 7. — A comparison of analytical (solid line) and numerical (filled circles) period spacings, as a function of Mcr/M«, where each has been normalized to 
the period spacing in the uncrystaUized case. In order to minimize mode trapping effects, the Schwarzscldld criterion has been used to compute the Brunt- Vaisala 
frequency in (a), whereas in (b) the modified Ledoux prescription has been used. The "kink" in (b) for 0.75 < Md/M* < 0.90 is due to the changing C/O profile in 
the core. 
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Fig. 8. — The same as Figure|jb, except that periods between 500 and 700 sec have been used to define the average period spacing from the pulsation calculations; 
we have picked this range of periods to mimic that observed in BPM 37093. For this case we see that mode trapping effects result in larger deviations from the 
asymptotic relation. 



K, artificially smooth by using the Schwarzschild criterion for 
the Brunt-Vaisala frequency, which essentially removes the 
bumps from the Brunt-Vaisala frequency and therefore mini- 
mizes mode trapping. The agreement between the two methods 
is extremely good. 

We now examine the more realistic case, where we include 
the modified Le doux criterion f or the Brunt-Vaisala frequency 
as described in [Brassard et al.| (1991). This plot is shown in 



Figure Mb. Although the overall shape of the plot has changed 
somewhat, the agreement between the asymptotic and numeri- 
cal results is still quite good. The observed "kink" for 0.75 < 
Mcr/M* < 0.90 is caused by the oxygen mass-fraction decreas- 
ing from 0.80 to 0.00 in this range. If we examine Figure ^ for 
the period formation region and we imagine moving the crys- 
tallization region to the right, we see that as we encounter the 
O/C transition zone, the rate of change of area under the curve 
doubles, so we would expect the slope of the curve in Figure^ 
to double as well, which is what we find. 

If we use a smaller range of periods to define the period spac- 
ing numerically, then we expect mode trapping effects to be 
amplified even further. This is illustrated in Figure ^, where we 
have used i = 2 periods in the range 500-700 sec to calculate 
a period spacing. Thus, if we have a complete set of observed 
^ = 2 periods in this range, we can typically expect "errors" of 
order ~5% in translating this to an asymptotic period spacing. 

An equivalent statement to the period spacing increasing 
with Mc,/M* is that the modes themselves are getting farther 
apart in period, so their periods must also be increasing. To il- 
lustrate this, we show how a spectrum of mode periods evolves 
continuously with Mcr/M*. Since mode identification between 
different models is not a simple matter, we have calculated the 
spectrum of modes on a fine enough mesh in M^^jM^ so that 
the period changes are small compared to the differences be- 



tween consecutive radial overtones. We then identify a given 
mode at one mesh point with the nearest mode in period of the 
neighboring mesh point. 

The result of this calculation for i = 2 periods is shown in 
Fi gure pL where the model considered is a 1.1 q with T^ff = 
11,800 K, MHe/M,^ = 10-^ and Mh/M,^ = 10"^ We have used 
the "hard-sphere" approximation for the solid/liquid bound- 
ary and the full Ledoux prescription for the Brunt-Vaisala fre- 
quency in calculating these periods. We have varied the param- 
eter Mcr/M* in increments of 0.01 from 0.00 to 0.99. We note 
that the periods either appear to be increasing or are relatively 
constant. In fact, even in regions in which the period of a given 
mode appears to be constant, its period is still slightly increas- 
ing with M^,lM^. 

Figure ^ represents the most detailed calculation to date 
showing how g-mode periods in white dwarf models evolve as 
a si ngle parameter is slowly varied (for a previous example, 
see Wood & Winget 1988). As such, it exhibits many inter- 



esting features. First, the "kink" in the periods in the range 
0.75 < Mcr/M* < 0.90 is again due to the crystallized region 
moving out into a region with a changing C/O profile; thus, 
this feature is merely a result of our assumed C/O profile. A 
more intrinsic feature of this plot are the "avoided crossings." 
While it is never possible numerically to establish with com- 
plete certainty whether or not true avoided crossings occur, the 
way in which the periods evolve as a function of M^/M^, is 
strongly r eminiscent of behavior found by A izenman, Smeyers, 
& Weigert (1977). For example, the two lowest period modes 
pictured in Figure ^ have what appears to be an avoided cross- 
ing at Mci/M,t = 0.58. To the left of this point, the lower period 
mode has more of its kinetic energy deep in the model near 
the solid/fluid interface, while to the right of this point it is the 
higher period mode which has its kinetic energy deeper. Thus, 
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the modes do switch character at thi s point, in the manner found 
by Aizenman, Smeyers, & Weigert (1977). 

Our general result that the ^-mode periods increase due to 
the presence of crystallization is not what was found by Hansen 
& Van Horn (1979), who reported that the ^-mode periods be- 
came shorter when the finite shear of the solid core was in- 
cluded. We believe that the resolution of this disagreement lies 
in a re-interpreta tion of their calculate d periods, not in the peri- 
ods themselves. Hansen & Van Horn (1979) calculated the pe- 
riods of ^ = 1 and 2 modes for £= 1,2, and 3, in both the fluid 
case and in the case of a 99.9% crystallized core. They found 
that in the crystallized case, the k = I periods had decreased by 
approximately a factor of two compared to the fluid case; for 
example, the £ = 1 period decreased from 193.8 sec to 99.8 sec. 
Our interpretation is that the 99.8 sec mode is actually a new 
mode, which would not exist if the core were not crystallized. 
Thus, the main effect of the solid core in their calculations was, 
in our view, to add an extra mode with a period below that of 
the previous ^ = 1 mode. To support this, we compare their ^ = 1 
periods in the fluid case with their k = 2 periods in the solid case. 
For £=l,2, and 3, we find that their periods now increase from 
193.8 to 193.9 sec, from 111.9 to 112.0 sec, and from 79.1 to 
79.2 sec, respectively. While these increases are small, they are 
consistent with what one might expect from a Tgff ^ 10,000 K 
Fe core white dwarf model which is strongly degenerate in its 
interior. In addition, the periods in the uncrystallized and the 
crystallized state are close enough to strengthen our conviction 
that this is actually the "correct" mode identification. 

Using our "global" code, we are numerically unable to treat 
models which are more than 97% crystallized. For 97% crystal- 
lized models, we do find evidence for low-period "interfacial" 
modes which do not exist in the uncrystallized case; interfa- 
cial modes such as these were foun d in neutron star models by 
McDermott, Van Horn, & Hansen! (1988). These m odes could 



be the new modes found by Hansen & Van Horn (1979). We 
caution, however, that we do not understand the properties of 
these modes, i.e., how they change period as the degree of crys- 
tallization changes and whether or not the standard definition 
of radial overtone number is still meaningful. We are therefore 
unable to extrapolate these res ults with confidence to the case 
of 99.9% crystallization which Hansen & Van Horn treated. 



7. (AP) AS A FUNCTION OF THE MODEL PARAMETERS 

In uncrystallized models, the period spacing is a function 
of many things, including the total stellar mass, the effective 
temperature, and the hydrogen layer mass. This is still true in 
the crystallized case, and we examine the effects which each 
has on (AP). The fiducial model against which we compare 
our calculations is a model with = 1 . IMq, Teff = 11, 800 K, 
Mh/M^ =10"^, and Mae/M^ = 10"\ Unless otherwise stated, 
all periods are calculated using the modified Ledoux prescrip- 
tion for the Brunt- Vaisala frequency. 

7.1. The Hydrogen Layer Mass, Mh 

For 0.6 Mq models, nuclear burning con siderations force 
Mjj/ Mi, to be smaller tha n a few times 10""* (Iben & Tutukov 



1984; [Iben & Macdonald) 1985). For models near 1.1 Mq, this 
translates into Mh/M* <^ 10"^ due to the higher gravities and 
pressures. We therefore examine models with Ma/Mi, between 
10"'° and 10"^ 

In Figure [To|a, we plot (AP) versus logMn/My, for differ- 
ent degrees of crystallization, as shown in the legend. For this 
model, we have used a C/O core and set MHe/M,^ = 10"^ and 



Teff = 11, 800 K, and we have calculated the Brunt-Vaisala fre- 
quency using the Schwarzschild criterion, so that we may min- 
imize mode trapping effects as much as possible. We see that 
the effect of increasing crystallization is to increase (AP) at all 
compositions. Similarly, the effect of decreasing logMe/M^ is 
also to increase (AP), for all degrees of crystallization. Thus, 
a change in one can mimic a change in the other. Figure [ic| b 
shows the more physical case where we have used the Ledoux 
prescription for calculating the Brunt-Vaisala frequency in this 
model. The same trends are still evident, but the period spacing 
itself has decreased by 3^ sec for all the models. This dif- 
ference is due to the nontrivial contribution of the composition 
transition zones. 

The horizontal dashed lines in Figures and are useful 
in demonstrating how observations could be used to constrain 
the parameter space. Here, we assume a "measured" value of 
(AP) ^ 17 sec. The dashed lines represent an uncertainty of 
5% in translating this "observed" (AP) to an asymptotic value, 
as is suggested by the deviations due to a finite sampling of the 
period range in Figure ^ For these calculations, we have calcu- 
lated the period spacing between consecutive 1 = 2 modes, all 
with m = 0. The exact same dependencies hold for the case of 
t=\ modes, if the mean period spacings are multiplied by a 
factor of \/3- 

From Figure [l0|b, we find the following constraints on our 
parameter space: -7 < logMn/Mj, < -5 and 0.00 < Mc,/M* < 
0.80. This is a fairly large range for each parameter, but they 
are now no longer independent. If we know one of them, 
then that can reduce the allowed range for the other. For in- 
stance, if the model is 50% crystallized, then we must have 
-6 < logMe/My, < -5. Additionally, from Figure [ic| a we see 
that there is no choice of parameters for which the period spac- 
ing matches the "observed" value. This demonstrates the large 
effect which the composition transition zones have on the aver- 
age period spacing. 

7.2. The Total Stellar Mass, 

We now consider models which differ only in mass from our 
fiducial model; all the other parameters are held fixed. In Fig- 
ure pl] a we plot the average period spacing for a set of = 
1.15 Mq models, again as a function of Mh, where we con- 
tinue to use the more physical Ledoux prescription for the the 
Brunt-Vaisala frequency. Since the more massive models are 
smaller in radius, they have a higher average density, and there- 
fore smaller periods and period spacings. For the less massive, 
1 .05 Mq models in Figure pT]b, we find the opposite is the case; 
these models are larger in radius and therefore have larger pe- 
riod spacings. 

7.3. The Effective Temperature, r^ff 

In Figure |2|we show how the mean period spacing for £ = 2 
modes varies as a function of the effective temperature of our 
fiducial models. The horizontal dotted lines again bracket an 
"observed" period spacing of 17 sec. We see that as the mod- 
els cool, the period spacing increases. This occurs because 
the models are becoming more degenerate. As the models ap- 
proach complete degeneracy, the Brunt-Vaisala frequency be- 
comes arbitrarily small, except in composition transition zones, 
so the periods and period spacings become large. 

7.4. Scaling Relations 
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The results of the previous sections may be used to obtain 
approximate scaHng relations for ( Af ) . Using models with Teff 
between 11,200 K and 12,800 K, between 1.05 Mq and 
1.15 Mq, and hydrogen layer masses with H = -logMn/M^t 
between 5 and 10 (all ranges inclusive), we obtain the follow- 
ing relation: 



( AP)fc2 = A/[ 1 + 0.54/(// - 5)]"-24m;' TTff , 



(6) 

whereA = 14.7 sec,M* =M*/(l.lM0),7;ff = W(12,000K), 
and / = (AP) / (Af )o is the ratio of the asymptotic period spac- 
ing at finite crystallization to that at zero crystallization, e.g., 
the solid line in Figure The bar on / indicates that we have 
chosen / for an "average" model, where by average we mean a 
model which has = 1 . \Mq , T^ff = 12, 000 K, Mh/M* = 10"^ 
andMHe/M,^= 10"^ 

Next, we examine the sensitivity of (AP) to small changes 
in these parameters. If we look at small variations around a 
model which is 50% crystallized and has a H layer thickness 
corresponding to // = 5, we find 



^ = 0.13...- + 0.15.//-1.70^-0.95^, 



(7) 



where we have defined Wcr = M^jMi,, and 6Y represents a 
small chang e in a given qu antity Y . From fits of spectra of 
BPM 37093, [Bergeron et al] (1995) find Teff = 1 1 , 740 ± 200 K, 
and Mi, = 1.09 ±0.05Mo. From this we see that the errors in 
the mass determination produce about 5 times the effect of the 
errors in the temperature determination. Thus, My, is the most 
important input parameter which the observations can provide. 
The quantities nid- and H, the crystallized mass fraction and the 
negative of the log of the hydrogen layer mass, respectively, 
are not observable quantities in the standard sense. They can 
only be determined from an asteroseismological analysis of a 
particular star, which leads us to the topic of the next section. 



8. MODE TRAPPING 

The traditional way to obtain information about the surface 
layer thicknesses of white dwarfs is to use mode trapping in- 
formation for individual modes, i.e., calculate AP^ = P^+i -Pk 
directly from the data set and match this to numerical calcula- 
tions. There is no reason why this will not work now, as long 
as we have enough well-identified consecutive overtones. 

In general, a transition zone may trap a mode in the region 
above it or below it. For a mode to be trapped in the outer hy- 
drogen layer, it needs to have a resonance with the He/H tran- 
sition region such that its vertical an d horizontal dis placements 
both have a node near this interface ( Brassard et al. 1992); this 
is the case which is traditionally referred to as mode trapping 
in the context of white dwarfs. If we imagine integrating this 
mode inward from the surface using the boundary conditions 
there, then we see that all this condition depends on is the mode 
frequency. Whether or not a frequency which would be trapped 
is indeed an allowable normal mode frequency does depend on 
the amount of crystallization in the core. From this, we see that 
it should be possible to disentangle the effects of crystallization 
and mode trapping^ 



More precisely. Brassard et al. (1992) find that the average 
period difference (AP), between successively trapped modes is 



(AP),= 



2^2 



Ndr/r 



rH 



(8) 



where rH is defined as the radius at the base of the hydrogen 
layer; the integral is therefore over the hydrogen surface layer 
only. We see that this does not depend on any of the proper- 
ties of crystallized region, but only on those of the hydrogen 
envelope. 
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In Figure |T3|, we plot the forward period difference, AP, = 
Pi+i - Pi, versus period, f,, for an equilibrium model with M^, = 
l.lMo,reff = 11,800 K, Mh/M* = 10"^ and Mee/M* = 10"l 
This shows how the mode trapping changes as the crystallized 
mass fraction is varied from 0.0 to 0.9 in increments of 0.1. In 
general, we see that the amplitude (strength) of the trapping de- 
creases with increasing crystallization. This is because as the 
degree of crystallization increases, all modes become more like 
envelope modes, which decreases the differences between the 
trapped and untrapped modes. 

Concerning the detailed structure of the mode trapping itself, 
the combined effect of the different transition zones makes it 
difficult to define a trapping cycle. Furthermore, we see that 
this structure changes significantly as the degree of crystalliza- 
tion is changed by only 10%. This suggests that we will need 
to examine the degree of crystallization in smaller increments. 

Figure [l^a is a more detailed version of Figure |l3| which 
shows how AP changes as the crystallized mass fraction is 
increased from 0.25 to 0.34 in increments of 0.01. First, we 
note that there are many trapping features which move uni- 
formly to the right as the the degree of crystallization is in- 
creased. For instance, there is a trapping feature with a pe- 
riod of ^ 580 sec at 25% crystallization which migrates to 
a period range of ~ 640 sec at 32% crystallization. Second, 
there are many features which remain relatively constant. The 
mode trapping structure in the range 420-500 sec is virtually 
unchanged, and the mode with a period of ^ 775 sec is also 
somewhat trapped in the majority of the panels. This ^ 775 sec 
mode has a period which does not evolve as rapidly as many 
of the other modes. Even so, its period changes by ~ 0.6 sec 
for every 5% change in the degree of crystallization. Since it 
is possible to measure pe riods to quite high accuracies of a few 
tenths of a second (e.g., Winget et al. 1991, 1994), we should 
in principle, using modes such as this as well as more sensitive 
modes, be able to derive quite accurate estimates of the crystal- 
lized mass fraction, if we are able to obtain a unique solution. 

Figure ^4|b shows a different range of crystallization, this 
time between 50% and 59%, for the same white dwarf model 
as was used for Figures 13 and [l4] a. Here, we see that there 
is again a trapping feature which migrates through the 600 sec 
region, as well as a mode with a period again at ^ 775 sec 
which tends to be trapped (at least in the upper 7 panels). The 
fact that some periods are strongly affected by changes in the 
crystallized mass fraction while others are not suggests that we 
may be able to disentangle the effects of different surface layer 
masses from those due to crystallization. 

9. OBJECTIVE FITTING PROCEDURES 

We need an automated procedure for searching parameter 
space, both to obtain more precise fits and to address the is- 
sue of uniqueness of fit. The sensitivity of the trapping features 
to the crystallized mass fraction is both a blessing and a bane: 
It is a blessing because this should allow us in principle to de- 
termine precise values for Mq /M*, and it is a bane because in 
practice it requires the computation of an enormous number of 
models on a fine grid in order to sample the parameter space 
adequately. 

We are currently exploring different methods which would 
address these issues. The first part of the problem is choosing a 
method such as "simulated annealing" or a "genetic algorithm" 
which can find global minima of multidimensional functions. 
The second part is automatically generating the equilibrium 
models with a given set of fit parameters, so that these mod- 



els can be examined pulsationally. Traditionally, the evolution 
of such models has been a "hands on" procedure, and this is 
true of our evolutionary models as well. 

Unfortunately, this problem is beyond the scope of this pa- 
per. Given the recent increase in both processor speeds and the 
degree to which problems are being parallelized, it should be 
possible to implement an objective fitting scheme which will al- 
low us to sample adequately the parameter space of the models. 
Such an approach is currently being developed at the University 
of Texas. This should, among other things, allow us to assess 
objectively the uniqueness of our asteroseismological fits. 

10. CONCLUSIONS 

In this paper we have calculated the effect which crystalliza- 
tion in the cores of our white dwarf models has on the frequency 
spectrum of pulsations. To a very high degree of accuracy, we 
find that the kinetic energy of the g-modes is excluded from 
the crystallized cores of our models. As the degree of crys- 
tallization is increased, the kinetic energy of these modes is 
pushed farther out from the central regions, and both the pe- 
riods and the mean period spacing (Af ) between consecutive 
radial overtones of the same £ increases. Using an "observed" 
value of (AP)f=2 =17 sec, we show how the range of possi- 
ble models can be constrained, and how mode trapping features 
may be used to obtain more precise information about these 
fits. Since some periods are strongly affected by changes in the 
crystallized mass fraction while others are not, we may be able 
to disentangle the effects of different surface layer masses from 
those due to crystallization. 

The introduction of a crystalline medium which is able to 
support shear does allow a new class of modes to exist, the tor- 
sional or f-modes. Since these modes have zero radial displace- 
ment, they should be unable to couple to the overlying fluid 
layers and should therefore remain unobservable. The p-modes 
have periods which are only a few percent different from their 
uncrystallized values. Since these modes are not observed to 
be excited, they are also not of interest in the context of BPM 
37093 and the other ZZ Ceti's. 

By investigating stars such as BPM 37093, asteroseismol- 
ogy may eventually be able to tell us whether crystallization 
occurs in the way we expect theoretically. Since crystallization 
and the effects of phase separation are the largest single sources 
of systematic uncertainties in the white dwarf luminosity func- 
tion, this would allow us to improve our estimates of the age 
of the Galactic disk as derived from the observed white dwarf 
luminosity function. Furthermore, since phase separation, if it 
occurs, affects the composition of the central regions of white 
dwarfs, this could systematically affect the observed charac- 
teristics of SNIa (e.g., total luminosity), which are believed to 
come from white dwarf progenitors. 

Finally, digital surveys now in progress promise to add con- 
siderably to the presently known number of cool white dwarfs. 
For instance, the Sloan Digital Sky Survey (3unn 1995) should 
increase this number by a factor of approximately 20, with the 
result that we may have 20 such stars with which to test the 
theory of crystallization. 
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Fig. 13. — APi = -P, (forward period difference) versus P, (period) for ^ = 2 modes. Each panel is labelled by the degree of crystallization assumed for the 
model, with the other model parameters being held constant. 
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Fig. 14. — The filled circles connected by hnes show the period spacing in the model versus the period, for degrees of crystaUization varying between (a) 25% 
and 34% and between (b) 50% and 59%. 
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APPENDIX 

SPHEROIDAL OSCILLATION EQUATIONS IN A CRYSTALLINE MEDIUM 

The Equations 

We define and ^/,, the radial and horizontal parts of the total displacement, in terms of the total vector displacement, i.e.. 



oO smO dm 



(Al) 



We take the following equations from Hansen & Van Horn (1979). The oscillation variables are 



Zl = — , 
r 



1 



Mo 



Z2 = — Xa + 2fi—— 



dr 



and 



6 

23 = — , 

r 



r r 



^4 = — 1 J 

Mo V 



where and are the radial and horizontal displacements, respectively, as defined in equation Al, r is the radius, A = Fi/^- ^ 
is the shear modulus, a = (''^^r) "^T'' with i= £(£+ 1). These variables are the same as those in [Hansen & Van Horn| (1979) 
except that we have divided zi and Z4 by /io = /i(r = 0), so that the equations are dimensionless.[] The fourth order system of equations 
(in the Cowling approximation) is then 



rZ3 
rz\ 



-( 1 + 2\5)zi + M0&2 + A&3 , 

1 i " 

— (-cT^pr^ - 4pgr+4nGp^r^ + 4/i/3(5)zi - 4/i&2 + — (pgr - 2p,l35)zz + &4 , 
Mo Mo 
MO 

-Z\ + Z4, 

M 



1 1(22 

(gpr-2pl35)z\-\5z2+ — ■^-pcr r +2p,5 



MO 



MO 



\{2l-\) + 2p{l-\) 



•Z3-3Z4, 



where the prime denotes J-, (5 = (A + 2^) \ (3 = 3A + 2/i, g is the acceleration due to gravity, and p is the density. 



Central Boundary Conditions 

Since the models we are considering are crystallized in the center, we need to obtain the boundary conditions in the center so 
that we may begin the outward integrations. If we assume that the solutions go like r' near the center, we find four solutions: 
s = £-2, £,-{£+ I), -(£ + 3). Only the first two solutions are regular at the origin, so they span the space of physical solutions. The 



{zi} = a 



; given by 


( ^ ) 




2(£-l) 


r'-^ + b 


l/£ 




V 2i£-l)/£ ) 





V 



ii+i)[\i+fj.{i-2)] 

2[A£(£+2)+/i(£2+2£-l)] 

(e+i)[\(f-e-3)+iJ.(e^-e-2)] 

2[xe(e+2)+tJ.(P+2e-\)] 
1 



where a and b are arbitrary coefficients and where p and A in the above formula ar e taken to have their central values. These two 
solutions for the eigenfunction near the center are equivalent to the relations given in Crossley (1975), if the Cowling approximation 
is used. 



' We note that there is a typographical error in the definition of za in Hansen & Van Horn ( 1 979) — an additional factor of 1 / r — but it does not propagate throughout 
the rest of their formulae. 
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The Solid/Fluid Interface 

In practice, we integrate each independent solution outward from the center With the exception of 23, the {z,} are continuous at 
the solid/fluid interface. Since Z4 = in the fluid, we choose the ratio of a and b such that za vanishes at this interface. This leaves 
only one overall normalization constant. Furthermore, y\ = z\ at the boundary. Since zi is also continuous, we have 

Z2 = Aa/ ^0 = AVgCyi -yi)/ Mo, 

where we have used the oscillation equations in the fluid to express a in terms of the Dziembowski variables {y,}. At the fluid/solid 
interface, if we solve for the {y,} in the fluid in terms of the {z,} in the solid then we find 

y\ = z\, 

Mo 

where Vg = gr/c^^, and A is now T\P since ji is zero in the fluid. Since we now have specified yi and y2 (up to an overall normal- 
ization constant which is present in the {z,}), we can now integrate the normal oscillation equations in the fluid (in the CowUng 
approximation) out to the photosphere of the model. 

The main difficulty i n applying this proce dure is that numerical noise can come to dominate the integrations in the crystalline 



core. The model which Hansen & Van Horq (1979) considered was a pure Fe core model near 10,000 K. As a result, the theory of 



crystallization suggested it should be about 99.9% crystallized by mass. The technique which we have used would probably not be 
viable for this case. The problem is that the two independent solutions, while quite different near the center, become almost linearly 
dependent farther out. Thus, we lose the ability to calculate the "difference" between the two solutions which is needed in order to 
set Za equal to zero at the solid/fluid interface. For our program, numerical noise dominates this process for g-modes in models which 
are more than 98% crystallized. 

In terms of the physics, however, we are somewhat over-dramatizing the situation, since nearly all of the pulsational results in 
this paper are based on the simple approximation that yi = at the solid/fluid boundary. From the self-consistent treatment, we have 
found this to be an extremely good approximation from 0% crystallization to 98% crystallization, and we have no reason to believe 
this situation will change at higher amounts of crystallization. Using this simplified treatment (yi = at the solid/fluid boundary), we 
are therefore able to treat accurately arbitrary degrees of crystallization. 

THE BOUNDARY CONDITIONS ON ^\R) AT THE CRYSTALLIZATION BOUNDARY 

The perturbations to the gravitational potential are, of course, generated by the perturbations in the density. They may therefore 
be written as 



<^>'(r) = G fdV'^, (Bl) 
J |r-r'| 



where we have assumed that the density vanishes at the surface. If we now write <l>'(r) = ^'{r)Y"\6,(j)) and p'jr' ) = p'{r')Y^ {9' 
and use the usual expansion of |r - r'| in surface harmonics, then we arrive at the result of jChristensen-Dalsgaard (1976): 



We will assume that there is no motion for r < r^, where is the radius of the crystaUization boundary. Thus, we have p'{r) = 
for r < rx. If we now take a derivative of <l>'(r) in the region < r < R, we find 

^ = ^ |-(^+l)r-(^-^> j^^'(r>'^-2j/ + ^/-'^%'(r')r'-^-'^r' + rpV)-rp'(r) 

where we have assumed that the density is zero at the outer boundary. Evaluating this at r = r^, and remembering that p'(r) = for 
r < Tx, we find that the first integral vanishes, which, along with the cancellation of the last two terms, leaves only 

£/$'(rx) 47rG 



dr le+i 



Ir'-' / p'ir'y'-'-'dr'. (B3) 



Finally, using equation B2 to evaluate <I>'(r) at r = r^, and combining it with equation B3, we obtain the final result 



'-^-U'^r^y ,B4, 

dr Tx 



This is the same boundary condition as is usually encountered in the uncrystallized case for r approaching zero (e.g., see 
1989); we see that it is unchanged by the presence of a rigid, crystallized core. 
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GLOBAL CONSERVATION OF MOMENTUM FOR t = 1 MODES 

As a note, we mention that the derivation by [Christensen-Dalsgaard (1976) that £ = 1 modes conserve the total momentum of the 
system and are therefore allowable pulsation modes is easily extended to the present case, in which we have a solid, completely rigid 
core, which is surrounded by a fluid envelope in which there are pulsations. 

In the present case, the center of mass of the fluid in the envelope is displaced by the pulsations. However, the pressure variations 
associated with these pulsations in the fluid exert a net force on the crystalline core, causing its center of mass to also move. 
Considered as a system, the core plus envelope conserves momentum, so that there is no net displacement of the center of mass, and 
i=\ oscillations are again allowed. 
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